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FOREWORD 


This  report  is  published  by  The  Aerospace  Corporation,  El  Segundo,  California, 
under  Air  Force  Contract  Nos.  F04701-70-C-0059  and  F04701-71-C-C172. 

This  report,  which  documents  research  carried  out  from  July  1970  through 
July  1971,  was  submitted  for  review  and  appx-oval  on  7  October  1971  to 
Capt.  James  L.  Warwick,  SAMSO/IND. 

The  Aerospace  Corporation  TRACE66  Trajectory  Analysis  and  Orbit  Deter¬ 
mination  Program  has  evolved  during  the  past  six  years  from  the  design  and 
implementation  efforts  of  many  individuals.  Recent  analysis  and  programming 
contributions  have  been  made  by  G.  Buechler,  W.  D.  Downs,  E.  H.  Fletcher, 
P.  T.  Gray,  and  A.  Z.  Rusick.  In  addition,  consultations  with  W.  T.  Kyner 
and  Lj.  Wong  have  led  to  many  significant  improvements  and  added  capabilities 
within  the  program.  The  TRACE66  documentation  series  is,  like  the  pro¬ 
gram's  growth  and  development,  a  joint  effort. 
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Director 
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NOTICE 


Some  of  the  previously  published  volumes  of  the  TRACE66  documentation 
series  have  been  published  by  The  Aerospace  Corporation  as  Technical 
Operating  Reports.  Volume  III:  Trajectory  Generation  Equations  and 
Methods  was  published  as  TOR-0066(9320)-2,  Vol  III;  Volume  V:  Differ¬ 
ential  Correction  Procedure  and  Techniques  as  TOR-0066(9320)-2,  Vol  V. 

Volume  VII:  Usage  Guide  was  published  as  TR-0O59(932O)-i,  Vol  VII,  and 
Volume  X:  Lunar  Gravity  Analysis  as  TR-0059(9320) -1,  Vol  X.  Future 
volumes  in  this  series  will  be  published  as  Technical  Reports. 

The  information  that  was  to  appear  in  Volume  VIII  is  included  in  Volume  VII; 
Volume  VIII  will  therefore  not  be  published. 


ABSTRACT 


The  TRACE66  Trajectory  Analysis  and  Orbit  Determination  Program  is  a 
general-purpose  orbital  analysis  program.  It  was  written  specifically  for 
the  CDC  6000  series  computers  to  assist  The  Aerospace  Corporation 
personnel  in  the  analysis  and  design  of  satellite  orbits  and  tracking  systems. 
Volume  I  is  a  reference  for  general  program  objectives;  it  describes  and 
summarizes  the  program.  A  comprehensive  description  of  its  areas  of 
application  is  presented.  An  overview  of  the  major  capabilities  of  the  program 
is  emphasized. 


The  TRACE66  documentation  series  is  summarized  as  follows: 


Volume  I: 
Volume  II: 
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Lunar  Gravity  Analysis 

LGA  Data  Processor 

Sequential  Least  Squares  &  Recursive  Filter  Procedures 
and  Techniques 


-v- 


CONTENTS 


NOTICE .  111 

ABSTRACT  . .  v 

1.  INTRODUCTION  AND  SUMMARY  .  1 

1.  1  Introduction  .  * 

1.2  Background  Information  .  2 

1.3  Applications  in  Orbital  Analysis  for  TRACE66  .  2 

1.3.1  Ephemeris  Generation  .  2 

1.3.2  Measurement  Data  Generation  .  3 

1.3.3  Measurement  Error  Modeling  .  3 

1.3.4  Orbit  Determination  .  3 

1.3.5  Covariance  Analysis  .  4 

1. 3.  6  Orbit  Determination  -  Postflight  Mode  .  5 

1.3.7  Improved  Modeling  .  5 

1.3.8  Multiple  Arc  Fitting  .  5 

1. 3.  9  Simultaneous  Vehicle  Mode  of  TRACE66  .  6 

1.3.10  Other  TRACE 6 6  Capabilities .  ? 

1.4  Structural  Characteristics  .  ' 

1.4.1  Functional  Characteristics  .  7 

1.4.2  Physical  Structure  . 

1.5  Documentation  Series  . 

2.  TRAJECTORY  GENERATION  OVERVIEW  .  15 

2.  1  Formulation  and  Numerical  Method  .  15 

2.2  Coordinate  and  Time-Keeping  Systems  .  17 

2.2.  1  Reference  Coordinate  Systems  . 1? 

2.2.2  Time-Keeping  Systems  . .  21 

2.3  Force  Model  . 23 

2.4  Trajectory  Partial  Derivatives  via  Variational 

Equations  .  25 


-vii- 


CONTENTS  (Continued) 


3.  EPHEMERIS  GENERATOR  OVERVIEW  .  27 

3.  1  Single  Vehicle  Ephemeris  Generation  .  27 

3.2  Simultaneous  Vehicle  Ephemeris  Generation  .  28 

4.  DATA  MEASUREMENT  GENERATOR  OVERVIEW .  31 

4.  1  Single  Vehicle  Data  Measurement  Generation  .  31 

4.2  Simultaneous  Vehicle  Data  Measurement 

Generation .  31 

5.  ORBIT  DETERMINATION  OVERVIEW  .  37 

5.  1  Batch  Differential  Correction  by  Weighted 

Least  Squares  .  38 

5.2  Sequential  Differential  Correction  by  Weighted 

Least  Squares  .  43 

5.3  Potter  Square -Root  Filter  .  44 

6.  COVARIANCE  ANALYSIS  OVERVIEW  .  45 

6.  1  Batch  Weighted  Least  Squares  Covariance 

Analysis  .  45 

6.  2  Sequential  Weighted  Least  Squares  Covariance 

Analysis  . .  47 

7  .  LUNAR  GRAVITY  FIELD  ANALYZER  OVERVIEW .  49 

REFERENCES  .  51 


-vm- 


FIGURES 


1.  Schematic  of  Major  Functions  .  8 

2.  TRACE66  Overlay  Structure  .  10 

3.  Accuracy  of  Numerical  Integration:  Analytic 

Solution,  Keplerian  Model  . . .  16 

4.  TRACE66  Coordinate  Systems  .  20 

5.  TRACE66  Time-Keeping  Systems  .  22 

6.  Force  Model  Representation  of  the  Equations 

of  Motion  .  24 

7.  Generation  of  Trajectory  Partials  . 26 

8.  Data  Measurement  Generator  Visibility 

Geometry  .  33 

9.  Schematic  of  the  Differential  Correction 

Procedure  .  40 

10.  Orbital  Statistics  via  Covariance  Analysis  .  46 

11.  Schematic  of  Sequential  Least  Squares  Covariance 

Analysis  Algorithm .  48 

TABLES 

1.  Representative  Numerical  Integration  Solution  Time 

Examples  via  CDC/6600  18 

2.  Single  Vehicle  Measurement  Types  .  34 

3.  Simultaneous  Vehicle  Measurement  Types  .  35 

4.  Orbit  Determination  Parameters  .  39 

-ix- 


1.  INTRODUCTION  AND  SUMMARY 


1.  1  INTRODUCTION 

This  document  is  the  first  volume  of  the  TRACE66  Trajectory  Analysis  and 
Orbit  Determination  Program  documentation  series.  It  is  directed  towards 
the  non-user  or  the  potential  user  interested  in  obtaining  an  overview  of  the 
capabilities  of  the  TRACE66  program. 

Section  2  discusses  the  characteristics  of  the  trajectory  integrator  with 
emphasis  on  the  following: 

•  Formulation  of  the  equations  of  motion  and  variational 
equations 

•  Numerical  technique  to  obtain  solutions 

•  Accuracy  of  technique 

•  Representative  solution  times 

•  Coordinate  and  time-keeping  systems 

•  Completeness  of  force  model 

Sections  3  and  4  discuss  the  ephemeris  and  data  measurement  generation 
functions  for  single  and  simultaneous  vehicle  modes. 

Section  5  presents  an  overview  of  the  orbit  determination  function  of 
TRACE66.  The  following  procedures  may  be  used: 

•  Batch  differential  correction  by  weighted  least  squares 

•  Sequential  least  squares 

•  Recursive  Potter  square -root  filter 

Section  6  describes  orbital  statistics  or  covariance  analysis  capabilities  of 
the  program,  and  Section  7  presents  an  overview  of  the  capabilities  of  the 
lunar  gravity  field  analyzer. 
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1.2 


BACKGROUND  INFORMATION 


TRACE66  is  a  computer  program  that  simulates  orbital  motion  and  tracking 
operations.  In  the  past,  orbit  determination  was  performed  at  The  Aerospace 
Corporation  by  the  TRACE-D  Orbit  Determination  Program  (Ref.  1),  which 
was  superseded  in  1967  by  TRACE66.  The  current  TRACE66  program  has 
capabilities  significantly  beyond  those  of  TRACE-D. 

1.3  APPLICATIONS  IN  ORBITAL  ANALYSIS  FOR  TRACE66 

TRACE66  is  The  Aerospace  Corporation's  trajectory  analysis  and  orbit 
determination  program;  its  applications  encompass  a  wide  range  of  problems 
in  orbital  mechanics.  In  general,  TRACE66  is  a  general-purpose  orbital 
analysis  program  used  to  assist  corporate  personnel  in  the  analysis  of 
tracking  operations  and  orbital  motion  of  artificial  satellites  about  the  earth, 
moon,  and  other  bodies  within  the  solar  system.  The  term  "orbit  deter¬ 
mination"  will  be  defined  and  discussed  later. 

A  hypothetical  situation  is  used  to  illustrate  some  typical  applications  of  the 
TRACE66  program.  Some  problems  not  addressed  by  TRACE66  are  also 
noted.  Suppose  you  are  the  manager  of  a  satellite  program  and  face  the 
problems  discussed  in  the  following  sections. 

1.3.1  Ephemeris  Generation 

Your  orbital  designers  have  somehow  (without  using  TRACE66,  which  does 
not  select  or  optimize  orbits)  picked  the  desired  orbit  for  the  first  flight. 

You  must  direct  your  recovery  ships  to  their  proper  stations. 

Ground  track  coordinates  are  some  of  the  outputs  of  ephemeris  generation, 
which  is  the  computation  of  position  and  velocity  components  and  related 
quantities  as  a  function  of  time.  Here,  the  required  initial  conditions  are 
those  of  the  nominal  orbit.  In  this  instance,  the  ephemeris  can  be  generated 
without  using  the  sophisticated  force  models  and  coordinate  systems  available 
in  TRACE66. 


1.  3.2 


Measurement  Data  Generation 


You  must  also  alert  the  tracking  network  and  provide  each  station  with  the 
predicted  rise  times  and  associated  azimuth  and  range. 

The  data  generation  capability  of  the  program  will  provide  the  desired  pre¬ 
dictions.  Here  again,  nominal  initial  conditions  and  simple  models  will 
suffice  for  the  accuracy  required.  As  a  little  plus  for  yourself,  get  the  rise- 
set-bar  chart  off  the  visibility  printer  plot.  It  will  show  at  a  glance  which 
stations  can  see  the  bird  at  any  given  time. 

1.3.3  Measurement  Error  Modeling 

For  their  dress  rehearsal,  the  tracking  analysts  need  some  realistic  simu¬ 
lated  tracking  data  for  exercising  the  operational  software  and  their  own 
analytic  techniques.  Another  data  generation  run  will  do  most  of  the  job. 
Biases  and  random  noise  should  be  added  to  the  generated  data,  and  perhaps 
the  environmental  model  (i.  e. ,  the  gravitational  model,  drag  coefficient, 
and  station  locations)  should  be  intentionally  falsified  to  provide  further 
realism.  The  simulated  data  will  be  in  TRACE66  format  (no  alternatives 
are  available),  but  it  may  be  placed  either  on  cards  or  on  magnetic  tape  for 
subsequent  conversion  to  another  format.  No  really  "wild"  points  will  be 
generated;  if  these  are  needed  to  test  editing  features  in  the  operational 
program,  manual  methods  will  be  required. 

1.3.4  Orbit  Determination 

The  launch  is  successful,  but  the  tracking  analysts  are  complaining.  Their 
orbit  determination  program  has,  according  to  plan,  adjusted  the  initial 
condition,  drag,  and  bias  parameters  so  that  the  computed  orbit  and  obser¬ 
vations  best  match  the  actual  data  received  from  the  tracking  sites,  but  too 
many  iterations  are  required,  and  many  of  the  final  discrepancies  (residuals) 
are  larger  than  the  radar  engineers  had  expected. 
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The  process  described  above,  the  estimation  of  unknown  parameters  from 
actual  observations,  is  the  orbit  determination  function  of  TRACE66.  The 
other  capabilities  of  TRACE66  are  provided  because  they  satisfy  naturally 
related  requirements  and,  in  many  instances,  use  the  same  basic  compu¬ 
tations.  TRACE66  is  well  suited  to  this  kind  of  problem.  Because  of  the 
great  flexibility  in  its  modeling,  it  can  simulate  the  computations  of  the 
operational  software.  Because  of  its  extensive  list  of  parameters  for 
differential  correction,  TRACE66  can  (in  the  hands  of  a  skilled  analyst) 
identify  previously  unsuspected  systematic  errors. 

The  first  step  is  to  place  the  tracking  data  on  magnetic  tape  for  efficient 
processing.  TRACE66  accepts  tracking  data  in  a  variety  of  formats,  so 
reformatting  is  probably  unnecessary.  One  iteration  of  orbit  determination 
provides  a  check  on  the  operational  program.  TRACE66  can  use  the  same 
model  and  biases  and  start  from  the  same  initial  conditions.  It  will 
(presumably)  generate  the  same  trajectory  and  derive  the  same  residuals  if 
the  prelaunch  checks  were  thorough.  The  printer  plot  of  residuals  may  dis¬ 
close  some  bad  observations  that  weren't  edited  by  the  operational  program. 
The  editor  in  TRACE66  isn't  perfect,  but  bad  observations  can  be  rejected 
on  subsequent  runs  even  if  they  are  buried  in  a  magnetic  tape. 

Suppose  the  problem  is  not  simply  bad  data,  and  the  printer  plot  of  residuals 
shows  an  obviously  systematic  pattern.  The  analyst  must  nominate  a  prob¬ 
able  cause,  although  this  may  be  difficult.  He  must  then  appropriately 
modify  the  model  (perhaps  also  specifying  further  parameters  for  differential 
correction)  and  try  again. 

1.3.5  Covariance  Analysis 

Since  reconstruction  and  prediction  are  not  going  well,  someone  reopens  the 
argument  for  rapid  deployment  of  the  spare  sensor.  Although  you  feel  that 
this  is  the  wrong  solution  to  the  problem,  or  perhaps  no  solution  at  all,  you 
promise  to  look  into  it. 


a  covariance  analysis  may  provide  some  objective  support  for  your  intuition. 

A  baseline  computation  will  estimate  the  potential  tracking  accuracy  of  the 
present  system.  Even  when  appropriate  uncertainties  in  station  locations 
are  factored  in,  the  estimate  will  be  optimistic  due  to  the  many  unmet 
assumptions  implicit  in  such  calculations.  The  analysis  will,  however,  show 
whether  significant  improvements  can  be  expected  if  another  sensor,  with 
its  attendant  random  errors  and  uncertainties  in  station  location,  time  bias, 
etc. ,  is  added  to  the  system. 

i .  3.  6  Orbit  Determination-Postflight  Mode 

The  early  tracking  problem  described  in  Sec.  1.  3.  4  has  been  solved  by  using 
TRACE66  in  its  postflight  mode  of  orbit  determination;  meanwhile,  the 
problem  has  vanished  from  the  real-time  operation. 

It  seems  that  the  cause  was  an  incompletely  spent  fuel  tank  that  was  supplying 
some  unexpected  thrust.  Good  fits  were  obtained  when  an  appropriate  thrust 
model  was  incorporated  into  the  trajectory,  and  the  rate  of  fuel  flow  was 
determined  by  differential  correction.  Residual  patterns  now  appear  random; 
more  important,  acquisition  predictions  are  being  met. 

1.3.7  Improved  Modeling 

Because  you  anticipate  similar  problems  on  subsequent  flights,  you  modify 
the  operational  program  to  include  a  thrust  model  and  the  appropriate  flow- 
rate  parameter. 

TRACE66  can  assist  in  checking  out  this  new  feature  by  generating  test 
tracking  data,  by  printing  a  trajectory  and  its  partial  derivatives  with  respect 
to  the  new  parameter  (for  comparison  purposes),  and  by  demonstrating 
graphically  (via  printer  plots)  the  correctness  of  the  partial  derivatives. 

1.3.8  Multiple  Arc  Fitting 

Although  operations  are  now  going  smoothly,  you  notice  that  the.  biases 
determined  for  one  of  your  stations  vary  erratically  from  pass  to  pass. 
Because  there  is  no  indication  of  faulty  hardware,  you  suspect  that  the 
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station  location  is  erroneous.  Large  quantities  of  data  are  available, 
including  some  taken  by  a  nearby  sensor  that  has  tracked  another  agency's 
satellites.  Some  careful  postflight  orbit  determination  is  required,  using 
a  variety  of  data  types  from  multiple  arcs  to  determine  common  parameters 
(in  this  case,  station  location). 

This  application  is  probably  TRACEob's  forte.  The  program  features  that 
contribute  to  this  analysis  are:  sophisticated  force  models  including  radiation 
pressure,  solar-lunar  perturbations,  and  a  wide  selection  of  models  of  the 
atmosphere  and  the  geopotential;  the  ability  to  handle  several  arcs  and  data 
types;  precise  coordinate  and  time-keeping  systems;  and  proper  use  of 
constraints  and  a  priori  estimates. 

The  proper  use  of  constraints  and  a  priori  estimates  merits  some  elabora¬ 
tion.  Presumably,  both  your  station  and  the  one  nearby  are  mislocated; 
their  separation  is  known  precisely,  but  their  absolute  locations  are  given 
only  approximately  (within  certain  standard  deviations).  The  constraint  on 
their  relative  locations  can  be  imposed  exactly,  and  the  previous  estimate 
can  be  weighted-in  with  the  current  solution.  Appendix  F  in  Vol.  VII  of 
this  documentation  series  (Ref.  2)  gives  some  other  examples  of  the  use  of 
the  constraint  capability  of  TRACE66. 

1.3.9  Simultaneous  Vehicle  Mode  of  TRACE66 

Proponents  of  a  new  navigational  satellite  program  claim  that  it  offers  more 
tracking  coverage  for  less  money.  By  operating  a  constellation  of  syn¬ 
chronous  satellites,  they  offer  to  furnish  ephemeris  data  on  their  birds  and 
a  black  box  to  be  carried  on  yours  that  will  return  both  direct  and  multipath 
ranges  between  their  vehicles  and  yours. 

The  simultaneous  vehicle  mode  of  TRACE66  can  be  applied  in  this  situation. 
Not  only  are  several  satellites  in  simultaneous  operation,  but  the  obser¬ 
vational  measurements  are  functions  of  the  positions  of  two  or  more. 
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Explicit  problems  similar  to  those  of  the  previous  sections  could  be  posed 
in  the  simultaneous  vehicle  environment,  but  for  this  introduction  it  is 
sufficient  to  note  that  all  four  of  the  typical  TRACE66  functions  (ephemeris 
generation,  data  generation,  orbit  determination,  and  covariance  analysis) 
can  be  exercised  in  this  mode  and  situation. 

1.3.10  Other  TRACE66  Capabilities 

This  hypothetical  situation  has  been  used  to  describe  briefly  a  few  of  the 
applications  of  TRACE66;  some  of  its  other  capabilities  are  listed  below: 

•  TRACE66  can  perform  its  four  basic  functions  using  lunar 
and  interplanetary  trajectories. 

•  In  both  earth-  and  moon-centered  applications,  the 
gravitational  attraction  can  be  computed  from  spherical- 
harmonic  or  point -mass  models. 

•  Orbit  determination  can  be  performed  using  a  batch 
least  squares,  sequential  least  squares,  or  recursive 
filter  algorichm. 

1.4  STRUCTURAL  CHARACTERISTICS 

1.4.1  Functional  Characteristics 

The  functional  structure  of  TRACE66  has  been  designed  and  developed  around 
its  four  major  applications,  which  are  orbit  determination,  ephemeris  gen¬ 
eration,  simulated  data  measurement  generation,  and  covariance  analysis. 
Figure  1  schematically  illustrates  the  functional  design  concept  of  the  pro¬ 
gram;  note  that  the  input  data  processor  and  the  trajectory  generator  are 
common  to  all  four  major  functions.  This  structure  allows  for  the  execution 
of  several  functions  sequentially  against  the  same  data  base.  To  illustrate 
the  potential  of  this  framework,  some  examples  and  usage  concepts  asso¬ 
ciated  with  the  functional  aspects  of  TRACE66  follow. 


In  a  typical  simple  case,  the  program  reads  all  input  data,  generates  the 
trajectory  file,  and  processes  it  for  »ne  desired  output.  This  may  be  a 
printed  ephemeris,  a  covariance  anaivsis,  or  simulated  tracking  data. 


More  complicated  cases  (or  a  series  of  cases)  can  also  be  executed.  The 
trajectory  file  may  be  simply  an  integration  from  given  initial  values,  or  it 
may  be  the  result  of  a  trajectory  reconstruction  from  observational  mea¬ 
surement  data;  in  fact,  both  kinds  of  trajectory  generations  may  appear  in  a 
single  case.  Several  processing  functions  may  be  executed  with  each  case. 
Finally,  cases  may  be  stacked  indefinitely  within  a  single  job  on  the 
computer. 

In  its  most  elementary  mode,  TRACE66  generates  or  processes  data  from 
a  single  vehicle.  The  observational  measurement  data  in  this  instance  per¬ 
tains  only  to  the  space  vehicle  and  its  relation  to  the  tracking  stations  or  to 
the  earth  itself. 

In  the  multiple  arc  mode,  TRACE66  can  generate  or  process  data  from 
several  space  vehicles,  but  each  data  item  is  associated  with  only  one  vehicle. 
The  vehicles  are  independent  and  need  not  be  in  orbit  simultaneously.  The 
data  may  be  associated  with  only  one  object,  yet  be  processed  separately  in 
different  arcs.  A  proper  solution  for  a  sensor  or  model  parameter  (for 
example,  a  station  location  or  a  gravitational  anomaly)  can  be  derived  from 
data  obtained  from  several  vehicles  during  different  time  periods.  Normally, 
the  multiple  arc  reconstruction  has  common  parameters,  either  naturally 
or  as  the  result  of  an  imposed  constraint;  otherwise,  the  reconstruction 
could  be  done  separately.  Note  that  the  multiple  arc  mode  is  a  simple  exten¬ 
sion  of  the  single  vehicle  mode  and  requires  no  special  identification  or 
handling. 

However,  when  observational  data  concern  the  position  or  velocity  of  two 
or  more  space  vehicles  orbiting  simultaneously,  the  program  is  run  in  its 
simultaneous  vehicle  mode.  The  program  handles  simultaneous  vehicle 
analyses  differently  from  the  way  it  handles  single  vehicle  or  multiple  arc 
problems.  Input  processing  and  numerical  integration  are  common  to  both 
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simultaneous  vehicle  and  single  vehicle  (multiple  arc)  modes.  However, 
the  four  functional  processors  indicated  in  Fig.  1  are  separate  and  inde¬ 
pendent,  depending  on  the  mode  of  program  operation. 

1.4.2  Physical  Structure 

TRACE66  is  written  entirely  in  CDC  6000  series  FORTRAN  and  runs  under 
the  SCOPE  operating  system.  In  its  standard  configuration,  the  program 
requires  approximately  34,  000  sixty -bit  words.  Various  numerical  limita¬ 
tions  (e.  g.  ,  350  spherical  harmonic  terms  o.r  100  ground  stations)  stated 
throughout  this  document  refer  to  this  standard  configuration.  The  majority 
of  these  restrictions  may  be  relaxed,  with  the  penalty  of  additional  core 
requirements.  On  the  other  hand,  the  program  may  be  run  with  full  sophis¬ 
tication,  but  reduced  capabilities,  in  as  few  as  approximately  25,  000 
sixty-bit  words. 

TRACE66  is  modular  in  design  and  employs  main,  primary,  and  secondary 
overlay  levels  to  perform  various  functions  (see  Fig.  2).  For  example,  the 
input  data  processor  and  the  integrator  use  all  three  levels;  the  ephemeris 
generator  resides  only  in  a  primary  overlay  level.  The  total  number  of 
modules  (overlays)  in  the  program  is  36;  it  is  not  expected  that  any  one 
application  would  use  all  the  overlays.  The  average  number  of  modules 
employed  for  particular  functions  is  tabulated  as  follows: 


FUNCTION 

NUMBER  OF  MODULES 

Orbit  Determination 

Covariance  Analysis 

Data  Generation 

Ephemeris  Generation 

10 

10 

8 

8 

TRACE66  is  a  file -oriented  system  and  uses  mass  storage  devices  (disk 
and/or  magnetic  tape)  extensively  for  communication  among  modules.  The 
program  has  significant  restart  capabilities  through  the  use  of  saved  files. 


In  addition,  varied  output  files  are  available  as  options;  they  may  then  be 
used  as  input  to  other  programs,  such  as  plotting  or  statistical  analysis 
programs. 

Much  of  the  flexibility  of  TRACE66  is  due  to  its  modular  structure  and  its 
file  orientation.  For  example,  significant  modifications  may  be  performed 
within  a  particular  overlay  with  little  or  no  effect  on  the  remainder  of  the 
program.  Furthermore,  one  may  add  new  functions  simply  as  a  separate 
overlay  and  thereby  utilize  the  applicable  capabilities  (e.  g. ,  the  input 
processor  or  the  numerical  integrator)  of  the  remainder  of  the  program. 
Because  of  this,  TRACE66  is  commonly  used  as  a  test  bed  for  new  algorithms 
and  techniques  and  as  a  simulator  of  other  orbit  determination  and  trajectory 
analysis  software. 

1.  5  DOCUMENTATION  SERIES 

The  TRACE66  documentation  series  is  summarized  as  follows; 

Volume  I;  General  Program  Objectives,  Description,  and  Summary  is 
directed  towards  the  non-user  or  potential  user  interested  in  obtaining  an 
overview  of  TRACE66  capabilities. 

Volume  II:  Coordinate  &  Time-Keeping  Systems  with  Associated  Trans¬ 
formations  is  a  technical  reference  for  the  coordinate  and  time-keeping 
systems  and  related  transformations  used  within  TRACE66. 

Volume  III:  Trajectory  Generation  Equations  St  Methods  (Ref.  3)  serves 
as  a  technical  reference  for  the  trajectory  generation  function  of  TRACE66. 

Volume  IV:  Measurement  Data  Generation  &  Observational  Measurement 
Partials  is  a  technical  reference  for  the  generation  functions  of  the  measure¬ 
ment  data  and  associated  observational  measurement  partial  derivatives  of 
TRACE66. 


Volume  V:  Differential  Correction  Procedure  and  Techniques  (Ref.  4) 
serves  as  a  technical  reference  for  the  batch  differential  correction  proce¬ 
dure  and  associated  techniques  used  within  TRACE66. 

Volume  VI:  Orbital  Statistics  Via  Covariance  Analysis  is  a  technical 
reference  for  the  orbital  statistics  generation  or  covariance  analysis  function 
of  TRACE66. 

Volume  VII:  Usage  Guide  (Ref.  2)  serves  as  a  reference  defining  all  input 
data  required  to  perform  any  of  the  TRACE66  functions.  Each  input  item 
is  defined,  and  all  basic  data  deck  structures  necessary  to  execute  TRACE66 
are  described.  Note,  however,  that  constant  changes  and  improvements  are 
being  made  to  the  program;  this  volume  may  not  include  a  description  of 
every  input  quantity  available  in  the  current  version  of  TRACE66. 

Volume  VIII:  Usage  Overview  will  not  be  published;  the  information  that  was 
to  appear  in  this  volume  has  been  included  in  Volume  VII. 

Volume  IX:  Detailed  Program  Structure  describes  the  program  structure  to 
the  subroutine  level. 

Volume  X;  Lunar  Gravity  Analysis  (Ref.  5)  serves  as  a  technical  reference 
for  the  Lunar  Gravity  Field  Analyzer  of  TRACE66. 

Volume  XI;  LGA  Data  Processor  serves  as  a  technical  reference  for  the 
LGA  data  processing  function  of  TRACE66. 

Volume  XII:  Sequential  Least  Squares  k  Recursive  Filter  Procedures  and 
Techniques  is  a  technical  reference  for  the  sequential  least  squares  and 
recursive  square-root  Potter  filter  procedures  and  associated  techniques 
used  within  TRACE66  to  perform  orbit  determination. 
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2.  TRAJECTORY  GENERATION  OVERVIEW 


The  trajectory  generation  function  of  TRACE66  is  common  to  all  program 
applications.  This  is  the  process  of  computing  the  position  and  velocity  of 
the  space  vehicle  as  a  function  of  time,  given  the  position  and  velocity  at  some 
initial  time  and  a  dynamic  model  that  accounts  for  the  accelerations  acting  on 
the  vehicle.  A  reference  frame  is  defined  within  which  the  components  of  the 
various  vector  quantities  are  expressed  and  the  accelerations  are  numerically 
integrated,  subject  to  the  given  initial  conditions.  This  section  presents  an 
overview'  of  the  formulation  and  numerical  integration  method  and  of  the 
reference  coordinate  systems  and  force  models  available  in  TRACE66. 

2.  1  FORMULATION  AND  NUMERICAL  METHOD 

TRACE66  numerically  integrates  the  total  vehicle  acceleration  vector  in  a 
Cartesian  coordinate  system  (Cowell  formulation  of  the  equations  of  motion). 

The  numerical  technique  used  to  perform  the  integration  is  a  predictor-corrector 
tenth-order,  Gauss- Jackson  differencing  scheme.  The  integrator  is  operated 
in  either  a  fixed-step  or  variable-step  mode  with  automatic  local  truncation 
error  controls.  A  fourth-order  Runge-Kutta  method  is  used  for  integrator 
starting  and  halving  procedures.  When  the  fixed-step  mode  is  used,  the  set 
of  equations  being  integrated  can  be  divided  into  two  groups.  The  second 
group  of  equations  may  be  integrated  at  a  step  size  that  is  a  multiple  (by  some 
positive  power  of  two)  of  the  step  size  used  for  the  first  group.  This  is  particu¬ 
larly  advantageous  for  an  orbit  determination  application,  in  which  variational 
equations  can  be  integrated  at  a  step  size  larger  than  that  of  the  equations  of 
motion  without  loss  of  accuracy. 

Time  (t)  is  the  normal  independent  variable  in  a  Cowell  formulation.  For 
most  problems,  precise  and  efficient  geodetic,  lunar,  and  interplanetary 
trajectories  are  obtained  using  the  variable-step  mode.  However,  in  the  case 
of  highly  eccentric  orbits,  it  has  been  demonstrated  that  regularization  of  the 
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independent  variable  leads  to  more  rapid  and  accurate  solutions  (see  Ref.  6). 
TRACE66  can  perform  its  numerical  integration  with  a  regularized  inde¬ 
pendent  variable  s  defined  such  that 


where  r  is  the  magnitude  of  the  vehicle  position  vector  at  time  t,  and  1  £  n  <  2. 

The  numerical  accuracy  of  the  TRACE66  integrator  is  documented  in  Ref.  3  and 
is  indicated  here  by  a  representative  example  of  RSS  (root-sum-square)  error 
obtained  by  differencing  the  numerical  solution  with  a  known  analytical  solution. 
The  results  are  shown  (see  Fig.  3)  as  a  time  history  of  RSS  error  on  a  log-log 
scale. 

Table  1  presents  representative  solution  times  for  the  numerical  integration. 
Note  that  the  eccentric  orbits,  cases  1  and  2,  have  comparable  results  for  the 
regularized  independent  variable  s. 

2.2  COORDINATE  AND  TIME-KEEPING  SYSTEMS 

2.  2.  1  Reference  Coordinate  Systems 

An  inertial  coordinate  frame  is  usually  chosen  as  the  reference  system  for  a 
problem  in  orbital  mechanics  so  that  the  equations  of  motion  are  in  a  simple 
form.  Regardless  of  the  system  chosen,  it  is  necessary  to  transform  to  other 
coordinate  systems;  for  example,  earth-fixed  coordinates  are  needed  to  evaluate 
the  geopotential  and  to  locate  trackers.  Therefore,  in  selecting  a  reference 
coordinate  system  and  associated  transformations,  a  compromise  is  made 
between  extreme  accuracy  and  ease  of  computation.  A  choice  completely  free 
of  known  errors  and  approximations  is  computationally  time-consuming  and  is 
unnecessary  for  many  applications.  TRACE66  allows  the  user  to  select  a 
reference  coordinate  frame  and  transformation  equations  according  to  the 
precision  required  for  the  particular  application.  For  example,  for  near-earth 
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Table  1.  Representative  Numerical  Integration  Solution 
Time  Examples  Via  CDC/6600 
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satellites  with  primarily  electronic  observations  recorded,  it  is  usually  not 
necessary  to  employ  rigorous  transformations  to  the  celestial  system.  (The 
celestial  system  is  used  for  the  tabulation  of  lunar,  solar,  and  planetary 
ephemerides;  star  catalogues;  and  as  a  reference  for  optical  measurements.  ) 
However,  for  precise  computations  over  extended  time  periods  and  for  the 
analysis  of  optical  measurements,  rigorous  transformations  must  be  used. 

All  TRACE66  coordinate  systems  originate  at  the  center  of  gravity  of  the  central 
body.  The  central  body  can  be  any  body  for  which  there  is  available  a  tabulated 
ephemeris  relative  to  the  earth.  The  fundamental  plane  and  the  principal  axis 
may  be  defined  as  parallel  to  any  of  the  following: 

•  True  earth  equator  of  date  (instant)  and  mean  equinox  of  midnight 
date  of  epoch 

•  Mean  earth  equator  and  mean  equinox  of  midnight  date  of  epoch 

•  Mean  earth  equator  and  mean  equinox  of  1950.  0  (celestial  system) 

•  Mean  earth  equator  and  mean  equinox  at  a  specified  base  date. 

The  first  of  these  systems  is  simplest  in  terms  of  the  transformation  matrices 
employed.  Precession  of  the  equator  and  equinox  from  midnight  date  of  epoch 
to  date  and  nutation  of  the  earth  equator  are  ignored  in  all  transformation 
equations.  Furthermore,  though  the  system  is  noninertial  (fixed  to  the  true 
earth  equator),  additional  terms  to  account  for  its  motion  are  not  included  in 
the  equations  of  motion  of  a  satellite. 

The  three  remaining  systems  are  inertial,  and  the  transformations  include  the 
effects  of  precession  and  nutation.  The  actual  transformation  matrices  are 
recomputed  at  a  user-supplied  time  interval.  Therefore,  it  is  still  possible  to 
trade  computation  time  for  accuracy  by  varying  the  length  of  the  re  computation 
interval.  For  example,  it  can  be  done  as  frequently  as  once  per  integration 
step  or  as  infrequently  as  once  during  the  numerical  integration  process.  In 
addition  to  considering  precession  and  nutation  effects,  TRACE66  can, 
optionally,  account  for  the  difference  between  the  instantaneous  pole  (or  spin 
axis)  of  the  central  body  and  a  point  attached  to  the  crust  by  computing  a  pole 
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TYPICAL  COORDINATE  TRANSFORMATIONS 
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Fig.  4.  TRACE66  Coordinate  Systems 


wander  matrix.  This  transformation  matrix  is  updated  at  the  same  time- 
frequency  as  the  precession  and  nutation  matrices.  Figure  4  indicates  the 
typical  coordinate  transformations  available  in  TRACE66. 

2.2.2  Time-Keeping  Systems 

The  foundations  of  dynamics  demand  that  "uniform"  time  exist  and  that  it 
correspond  identically  to  the  time  variable  in  dynamical  equations.  Any 
physical  measure  of  time  is  established  by  definition;  it  usually  represents  an 
attempt  to  approximate  uniform  time  more  closely  than  did  its  predecessor, 
i.  e„ ,  to  reduce  discrepancies  between  observations  and  dynamic  theories. 

Though  each  of  the  definitions  is  nonuniform  to  some  degree,  many  have  areas 
of  continuing  application.  Four  time-keeping  systems  may  be  used  simulta.n- 
eously  in  TRACE66  to  solve  a  problem  in  orbital  dynamics.  These  systems  are: 

•  IT:  Uniform  integration  time  is  used  as  the  independent  variable 

in  the  equations  of  motion.  A  typical  example  is  IT  =  A1 
(atomic  time). 

•  ET:  Ephemeris  time  is  needed  to  relate  to  the  planetary 

ephemerides. 

•  UT1:  Universal  time  is  based  on  the  rotation  rate  of  the  earth  and 

is  required  to  compute  the  sidereal  angle. 

•  OT:  Observation  time  is  used  as  the  time  tag  identifying  the 

instant  at  which  a  measurement  is  recorded  at  a  tracking 
station.  A  typical  example  is  OT  =  UTC  (broadcast  time). 

It  is  assumed  in  TRACE66  that  IT  and  A1  are  related  to  ET  by  constants  and  that 
UT1  and  OT  are  both  related  to  A1  by  quadratic  polynomials,  with  time  in 
seconds  from  a  fixed  epoch  as  the  independent  variable.  Complete  flexibility 
is  afforded  the  user  by  varying  input.  For  example,  measurements  may  be 
time -tagged  in  any  time  system  by  using  the  polynomial  relating  OT  to  Al. 

For  a  short-duration  run,  it  can  be  assumed  that  UT1  is  uniform;  UT1  can 
then  be  used  as  the  independent  variable  for  the  equations  of  motion  if  the 
relationships  between  IT  and  ET,  ET  and  Al,  and  Al  and  UT  are  set  to  the 
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Fig.  5.  TRACE66  Time-Keeping  Systems 


appropriate  constants.  Figure  5  shows  a  typical  set  of  relationships  among 
time-keeping  systems.  Reference  7  contains  a  more  complete  description  of 
this  subject. 

2.  3  FORCE  MODEL 

The  f^rce  model  representation  within  TRACE66  can  range  from  a  simple 
central  gravitational  term  to  a  highly  sophisticated  model.  For  example,  the 
motion  of  the  vehicle  can  be  considered  to  be  influenced  by  the  following 
perturbative  effects: 

•  Gravitational  potential  of  the  central  body  expressed  as  spherical 

harmonics,  using  up  to  350  terms  of  C  ,S 

6  r  nm  nm 

•  Up  to  50  point  masses  (p.j,  pj,,  \j)  imbedded  in  the  central  body 

•  Gravitational  attraction  from  other  bodies  in  the  solar  system 

•  Atmospheric  density 

•  Ballistic  coefficient(s): 

*  Constant 

*  Time -segmented 

*  Represented  as  a  polynomial  in  time 

*  Represented  as  a  tabular  function  of  altitude,  time, 
Mach  number,  or  angle -of-attack 

•  Density  models: 

*  ARDC  1959 

*  U.S.  Standard  1962 

*  Lockheed-Jacchia 

*  Exponential 

*  LMSC  1967 

*  Jacchia  1964  (log  p,  Walker  analytic,  or  Walker 
modified  by  Bruce) 

*  AFCRL  (Champion  1968) 

*  NWL 

•  Seven  finite  thrust  models  with  achieved  velocity  increment  cutoff 
option 
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Fig.  6.  Force  Model  Representation  of  the  Equations  of  Motion 


•  Solar  radiation  pressure  constant  with  a  shadow-dependent  scale 
factor 

•  Orbit  adjusts  at  specified  times  by  instantaneous  change  in  velocity 
The  total  acceleration  a.cting  on  the  vehicle  is  given  symbolically  a.s 


r 


X  h 

1=0 


This  representation  includes  central  gravitation  and  all  the  perturbative  effects 
defined  in  Fig.  6.  Each  of  the  acceleration  components  ih  is  evaluated  in  its 
appropriate  reference  frame,  rotated  if  .necessary,  and  accumulated  in  a 
body-centered  rectangular  coordinate  system. 

2.  4  TRAJECTORY  PARTIAL  DERIVATIVES  VIA  VARIATIONAL 

EQUATIONS 

The  solution  of  an  orbit  determination  or  covariance  analysis  problem  requires 
the  ability  to  predict  the  effects  of  perturbations  of  various  parameters  on  the 
orbit.  TRACE66  estimates  these  effects  by  numerically  integrating  the  partial 
derivatives  of  the  total  acceleration  with  respect  to  the  parameter  of  interest. 
These  partial  derivatives  are  called  variational  equations,  and  the  resulting 
integrals  (solutions  to  the  variational  equations)  are  the  trajectory  partial 
derivatives.  Functionally,  the  total  vehicle  acceleration  vector  can  be 
expressed  as 


£  =  £.(£>  £>  P) 

where  p  is  an  explicit  or  implicit  equation-of-motion  parameter.  Figure  7 
symbolically  illustrates  (employing  the  above  functional  notation)  the  formation 
of  a  variational  equation  and  its  solutions. 
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p  IS  AN  EQUATION  OF  MOTION  PARAMETER 


3.  EPHEMERIS  GENERATOR  OVERVIEW 


The  ephemeris  generation  function  of  TRACE66  produces  vehicle  trajectory- 
information  at  prespecified  or  requested  event  times.  The  integrator  gener¬ 
ates  the  vehicle  state  S  (r,  £,  ii)  at  integration  times;  to  obtain  the  vehicle 
state  S(r_,  r ;  t)  at  a  particular  time,  Hermitian  interpolation  formulas  are  used 
(see  Ref.  8).  This  print  time  could  result  from  prespecified  print  times 
such  as: 


•  Epoch  event 

•  Print  time  vector 

•  Stop  event 


or  from  the  occurrence  of  trajectory  events  such  as: 

•  Instantaneous  orbit  adjusts 

•  Ballistic  coefficient  changes 

•  Start  and  stop  of  thrust  periods 

•  Observational  measurement  times 

•  Crash  altitude  event 

•  Particular  latitudes,  longitudes,  or  altitudes 

•  Ascending  and  descending  nodes 

•  Apsis  points 

•  Eclipsing  (entry /exit)  points 

•  Points  of  closest  approach 

•  Local  moon  or  midnight 


'v 


f  specified 


j 


-\ 


J-  detected 


3.  1  SINGLE  VEHICLE  EPHEMERIS  GENERATION 


The  types  of  trajectory  information  available  in  TRACE66  for  a  single  vehicle  are 
described  in  this  section. 
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The  vehicle  state  vector  is  characterized  as: 


Bodycentric 

x,  y,  z,  x,  y,  z 

Body-fixed 

.  ^  AAA 

A  A  A  *  •  » 

x,  y,  z,  x,  y,  z 

Spherical 

a ,  6,  (3,  Az,  R,  V 

Orbit-plane 

R,  T,  G,  R,  T,  C 

Classical 

a,  e,  i,  £2,  co,  t 

Related  quantities  are  represented  as: 


Latitude,  longitude,  and  altitude 

$,  X,  h 

Subvehicle  latitude  and  geomagnetic 
latitude  and  longitude 

**  ®g>  h 

Rev  number  and  apocenter  and 
pericenter  altitudes 

Rev.  ha,  hp 

Empirical  nodal  regression, 
period,  and  period  decay 

6'  P'  p' 

’  £2’  £2 

Secular  rates  of  node  and  pericenter 

£2,  cj 

Eccentric,  mean,  and  true  anomalies 

E,  M,  v 

Keplerian,  anomalistic,  and  nodal 
periods 

PK’  PA’  PN 

Right  ascension  and  declination  of 
planetary  bodies  (e.  g. ,  the  sun  or 
moon) 

<*',  6' 

Relative  position  and  velocity  vectors 
to  planetary  bodies  and/or  to  another 
vehicle  (e.  g. ,  trajectory  differences) 

•<JI 

<fl 

SIMULTANEOUS  VEHICLE  EPHEMERIS  GENERATION 

Frequently,  it  is  necessary  to  analyze  trajectory  information  for  a  vehicle 
system  that  consists  of  more  than  one  vehicle.  A  convenient  way  of  presenting 
this  information  is  to  time-collate  the  vehicle  ephemeris  output  at  prespecified 
print  times  that  result  from  the  following: 

•  Epoch  event 

•  Print  time  vector 

•  Crash  altitude  event 

•  Stop  event 
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This  time -collated  output  is  intermixed  with  individual  vehicle  ephemeris 
output  (see  Sec.  3.  1)  at  specified  trajectory  events  (thrust  start/ stop,  orbit 
adjusts,  ballistic  coefficient  change,  etc.  )  and  at  detected  trajectory  events 
(nodes,  apsis  points,  etc.  ). 

Optional  information  available  in  TRACE 66  during  the  use  of  this  function 
consists  of  the  following: 

•  Trajectory  differences  between  reference  vehicles  (maximum  of  five) 
and  the  remaining  vehicles  of  the  system. 

•  Vehicle -to -vehicle  visibility  information 

•  Visibility  status  array 

•  Range  and  range -rate  measurements 

•  Station-to-vehicle  visibility  information 

•  Visibility  status  array 

•  Range,  elevation,  and  azimuth  measurements 

•  Vehicle -to -vehicle  encounter  detection 

•  Time  and  distance  of  the  point  of  closest  approach  between 
two  vehicles 

The  basic  restrictions  and  limitations  associated  with  the  generation  of 
simultaneous  vehicle  ephemeris  information  are: 

•  Maximum  of  20  vehicles. 

•  Maximum  of  5  reference  vehicles  allowed  within  a  system  of  n 
vehicles  when  orbit  differences  are  computed. 

•  Epoch-stop  time  interval  for  the  first  vehicle  must  span  (or  have 
coincident  end  points)  with  all  other  vehicle  epoch-stop  time  intervals. 

•  All  trajectories  must  be  generated  by  an  earth-centered-inertial 
(ECI)  integrator. 
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4.  DATA  MEASUREMENT  GENERATOR  OVERVIEW 


The  data  measurement  generation  function  of  TRACE66  produces  simulated 
observational  measurements  and  visibility  information.  Hermitian  interpol¬ 
ation  formulas  are  used  to  obtain  the  vehicle  state  S  (r,  ty  t)  at  a  particular 
measurement  time. 

4.  1  SINGLE  VEHICLE  DATA  MEASUREMENT  GENERATION 

The  single  vehicle  measurement  data  generator  is  characterized  by  the  following: 

•  38  single  vehicle  measurement  types  (see  Table  2) 

•  Maximum  of  100  stations 

•  Vehicle-to-station  rise,  n.aximum  elevation,  and  set  detection 

•  Vehicle -planetary  body  (or  earth-fixed  station)  occultation  detection 

•  Station  visibility  constraints  on  range,  elevation,  and  azimuth 

•  Simulated  measurement  noise  and  bias 

•  Simulated  refraction  and  light-time  corrections 

Figure  8  shows  a  visibility  example  of  station -to -vehicle  rise /set  occurrence 
based  on  minimal  station  elevation  constraint  and  vector  geometry  related  to 
the  body,  stations,  and  vehicle. 

4.2  SIMULTANEOUS  VEHICLE  DATA  MEASUREMENT  GENERATION 

The  simultaneous  vehicle  measurement  data  generator  is  characterized  by 
the  following: 

•  19  simultaneous  vehicle  measurement  types  (see  Table  3) 

•  Maximum  of  100  stations 

«  Station  visibility  constraints  (e.  g. ,  minimum  and  maximum  range, 
azimuth,  and  elevation) 

•  Simulated  measurement  noise,  white  or  colored 

•  Simulated  station -dependent  and  vehicle -dependent  measurement 
biases 
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Preceding  page  blank 


•  Simulated  station-clock  error  models  with  both  deterministic 
and  statistical  terms 

•  Simulated  tropospheric  and  ionospheric  refraction  corrections  for 
some  measurement  types 

•  Simulated  light-time  corrections  for  some  measurement  types 

•  Station -to -vehicle  and  vehicle -to -vehicle  visibility  status  as  a 
function  of  time  is  available  in  matrix  and  printer  plot  formats 

•  Mutual  visibility  status  involving  sets  of  two  prespecified  station- 
satellite  combinations. 
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Table  2.  Single  Vehicle  Measurement  Types 


Symbol 

Definition 

Symbol 

Definition 

R,  R,  R 

Slant  Range,  Rate, 

0  ) 

Acceleration 

o 

Aspect  Angles 

Az,  Az,  Az 

Azimuth  Angle,  Rate, 

J 

Acceleration 

SA 

Signal  Attenuation 

Ei,  Ei,  Ei 

Elevation  Angle,  Rate, 

x  \ 

Acceleration 

1 

v  i 

Earth-Fixed  Rectan- 

P 

gular  Coordinates 

6 

2  ) 

"O 

>- 

Interferometer  Data 

aT  ) 

Topocentric  Right 

ir 

S 

Ascension  and 

Q 

6t  j 

Declination 

* 

Q'r 

Geodetic  Vehicle 

a  ! 

Geocentric  Right  Ascen- 

Latitude 

sion  and  Declination 

X 

Vehicle  Longitude 

HA 

Topocentric  Hour  Angle 

SR 

Surface  Range  from 

u 

Vehicle- Centered 

Argument  of  Latitude 

H 

Vehicle  Height 

V 

Vehicle  Cross-Plane 

D 

Doppler  Rate 

Angle 

f2 

Two-Way  Doppler 

:*  i 

Antenna  x  and  y  Angles 

fT 

Tranet  Doppler 

y  / 

• 

rsgls 

SGLS  Range  Rate 

A 

Geoceiver  Range 

Differences 

LA 

Look-Angle 

K 

Angie  between  Station- 

Line-of-Sight  and 

Geocentric  Vehicle 

Radius  Vector 

Table  3.  Simultaneous  Vehicle  Measurement  Types 


Symbol 

Definition 

R,  R 

Station-to- Vehicle  Slant  Range  and  Range  Rate 

rsgls 

SGLS  Range  Rate  (Range  Differences) 

Az 

Station-to- Vehicle  Azimuth  Angle 

El 

Station-to-Vehicle  Elevation  Angle 

Ax  | 

Station-to-Vehicle  Antenna  x  and  y  Angles 

Ay  ) 

V2,  V2 

Vehicle-to-Vehicle  Range  and  Range  Rate 

S2,  S2 

Station-to- Vehicle-to -Vehicle  Range  and  Range  Rate  Sums 

V3,  V3 

Vehicle-to-Vehicle-to-Vehicle  Range  and  Range  Rate  Sums 

S3,  S3 

Station-to-Vehicle-to-Vehicle-to-Vehicle  Range  and 

Range  Rate  Sums 

TOA 

Time  of  Arrival 

TDOA 

Time  Difference  of  Arrival 

3WR 

Three-Way  Range 

MP 

Multipath 

5.  ORBIT  DETERMINATION  OVERVIEW 


The  basic  orbit  determination  problem  is  to  estimate  values  for  a  set  of 
observational  and  model  parameters  that,  in  some  way,  best  represent  a 
given  set  of  measurements.  The  choice  of  parameters  is  dictated  by  the 
purpose  and  objectives  of  the  user  and,  to  some  extent,  by  the  nature  of  the 
available  measurements.  The  two  criteria  used  in  TRACE66  for  determining 
parameter  values  that  represent  observables  are: 

•  Minimization  of  a  prespecified  cost  function 

•  Production  of  a  minimum-variance  (optimal)  estimate 

The  first  criterion  is  empirical,  whereas  the  minimum-variance  approach  is 
formulated  in  statistical  terms.  Additional  empirical  devices  are  available 
in  the  program  to  assist  the  user  in  achieving  the  desired  minimum  (see  Sec. 

5.  1).  In  most  instances,  the  empirical  approach  is  used  in  solving  an  orbit 
determination  problem.  However,  if  certain  assumptions  are  made  (par¬ 
ticularly  with  regard  to  the  statistical  significance  of  contributors  to  the  cost 
function),  it  can  be  shown  that  the  minimization  of  a  cost  function  produces  a 
minimum-variance  estimate. 

The  second  criterion  used  for  estimation  allows  statistical  considerations  to 
be  imposed  on  the  solution  with  either  a  recursive  or  nonrecursive  algorithm. 

In  TRACE66,  three  different  computational  algorithms  are  available  to  process 
a  given  set  of  measurements  (Volumes  V  and  XII  of  this  documentation  series 
are  the  technical  manuals  for  these  algorithms).  These  three  orbit  determi¬ 
nation  algorithms  are  listed  as  follows: 

•  Batch  differential  correction  by  weighted  least  squares 

•  Sequential  differential  correction  by  weighted  least  squares 

•  Potter  square-root  filter 
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The  first  two  algorithms  are  nonrecursive;  that  is,  they  update  parameter 
values  at  one  time  point  based  on  a  set  of  data  that  may  be  collected  over  a 
relatively  long  time-arc.  One  obtains  an  estimate  at  any  other  time  point  by 
integrating  the  trajectory,  using  the  updated  parameter  values  as  initial 
conditions.  The  third  algorithm  is  an  example  of  recursive  estimation,  or 
filtering,  in  which  a  set  of  measurements  is  processed  successively  in  its 


natural  time  order,  and  an  estimate  is  produced  at  each  measurement  time. 

All  three  estimation  algorithms  are  described  in  the  following  sections,  with 
emphasis  on  the  following: 

•  Basic  equations  (cost  functions)  and  methods 

•  Convergence  characteristics 

•  Parameter  constraints 

•  Statistical  considerations. 

The  limitations  in  the  application  of  these  algorithms  are  listed  below: 

•  Single-  or  multi-vehicle  fit  with  no  limit  on  the  number  of 
vehicles  for  single-vehicle  measurement  types.  There  is  a 
limit  of  20  vehicles  when  simultaneous  vehicle  types  are 
considered. 

•  Total  of  100  parameters 

•  Maximum  of  60  equation- of- motion  parameters  for  any  vehicle 

The  orbit  determination  parameters  available  with  these  algorithms  are 
summarized  in  Table  4. 

5.  1  BATCH  DIFFERENTIAL  CORRECTION  BY 

WEIGHTED  LEAST  SQUARES 

The  criterion  used  in  this  algorithm  minimizes  the  differences  (residuals) 
between  the  actual  measurements  and  the  corresponding  values  computed  from 
the  mathematical  model  in  a  generalized  least  squares  sense  By  using  a 
weighting  matrix,  one  may  assign  the  proper  relative  importance  to  measure¬ 
ments  of  various  types  and  qualities.  Measurements  are  usually  related  to 
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equation- of- motion  and  observational  parameters  by  complicated  and  highly 
nonlinear  functions.  To  make  the  minimization  computationally  tractable, 
one  usually  linearizes  the  problem  and  uses  an  iterative  technique. 

An  initial  value  for  each  parameter  is  required  to  start  the  procedure.  Each 
iteration  solves  for  a  set  of  corrections  to  the  currently  estimated  values. 
The  computational  steps  involved  in  batch  weighted  least  squares  are  illus¬ 
trated  in  Fig.  9.  Note  that  in  this  context,  the  word  "batch"  implies  the 
processing  of  all  available  measurements  on  each  iteration  of  the  differential 
correction  procedure. 

The  basic  equation  of  weighted  least  squares  is  written: 

[ATWA]A.P  =  A1  WO 

J  me 
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where 


A  =  the  matrix  of  partial  derivatives  of  measurement  with  respect 
to  the  parameters 

W  =  the  weighting  matrix,  which  is  assumed  to  be  diagonal  or  block- 
diagonal  with  block  dimensions  <60 

O  =  the  vector  of  measurement  residuals 
me 

AP  =  the  vector  of  parameter  corrections  to  be  determined 

This  matrix  equation  is  referred  to  as  the  normal  equation,  and  the  matrix 
T 

[AaWA]  is  called  the  normal  matrix. 

When  weighted  least  squares  is  applied  to  an  orbit  determination  problem,  it 
is  not  uncommon  for  the  iterative  process  to  exhibit  poor  convergence,  or 
even  divergent  characteristics.  The  following  conditions  could  lead  to  such 
circumstances:  inadequacies  in  the  observational  model  and/or  a  poor  initial 
approximation  of  the  parameter  values.  In  such  situations,  the  analyst  may 
aid  in  the  eventual  convergence  of  the  differential  correction  process  by  im¬ 
posing  a  side  condition  on  the  solution  that  bounds  the  magnitude  of  the  correc¬ 
tion  vector  on  any  iteration.  The  bounds,  which  are  specified  individually 
for  each  parameter,  are  adjusted  dynamically  as  the  successive  iterations 
converge  or  diverge. 

A  related  aid  to  convergence  is  the  addition  of  a  constant  conditioning  matrix 
to  the  least  squares  normal  matrix.  A  good  choice  for  the  conditioning  matrix 
is  an  a  priori  parameter  variance/ covariance  matrix.  Note  that  both  of  these 
convergence  aids  are  strictly  empirical  and  have  no  statistical  ramifications. 
Since  the  cost  function  being  minimized  is  unchanged,  the  final  parameter 
values  are  unrestricted  and  will  coincide  with  the  ordinary  weighted  least 
squares  solution. 

If  one  assumes  that  the  mathematical  model  is  exact,  that  the  observational 
errors  are  random  with  mean  zero  and  variance/covariance  matrix  W  ,  and 
that  the  observations  are  linear  functions  of  the  parameters,  it  can  be  shown 
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that  the  weighted  least  squares  solution  is  a  minimum-variance  solution  (see 

T  -1 

Ref.  9).  Furthermore,  the  inverse  normal  matrix  (A  WA)  equals  the 
parameter  error  variance/ covariance  matrix.  This  matrix  is  computed;  the 
resulting  TRACE66  output  serves  as  an  estimate  of  the  reliability  of  the  orbit 
determination  results. 

When  the  initial  parameter  values  are  not  merely  guesses,  it  may  be  desir¬ 
able  to  limit  the  parameter  corrections  according  to  the  confidence  level 
associated  with  the  initial  estimate.  The  reason  for  limiting  these  corrections 
may  stem  entirely  from  empirical  considerations,  based  on  confidence  in 
past  values  and  the  continuity  of  the  process.  Or  the  reason  might  stem  from 
statistical  considerations;  in  this  case  the  objective  is  to  statistically  update 
parameter  estimates  according  to  physical  measurements  and  their  associated 
statistics.  With  either  interpretation,  a  slight  modification  to  the  basic  equa¬ 
tion  of  weighted  least  squares  is  required. 

When  one  views  orbit  determination  by  the  batch  weighted  least  squares 
method  as  an  empirical  process,  one  may  use  either  of  two  approaches  to 
derive  the  applicable  equations.  The  first  approach  is  to  augment  the  mea¬ 
surement  vector  and  its  weighting  matrix  with  the  initial  estimate  and  an 
appropriate  set  of  weights.  The  second  approach  is  to  redefine  the  cost 
function  as  the  weighted  sum  of  squares  of  the  measurement  residuals  plus 
the  weighted  sum  of  squares  of  the  total  parameter  corrections.  Either  way, 
one  obtains  the  same  equations. 

On  the  other  hand,  if  one  regards  the  orbit  determination  process  as  statis¬ 
tical,  the  problem  is  that  of  optimal  combination  of  statistical  estimates 
(i.  e.  ,  an  a  priori  estimate  and  its  covariance  matrix  are  to  be  combined  with 
the  estimate  and  covariance  matrix  produced  by  the  previously  described 
minimum- variance  algorithm).  If  optimal  combination  is  defined  as  a  linear, 
unbiased,  minimum- variance  combination,  the  resulting  equations  correspond 
to  those  of  the  empirical  algorithm,  with  the  inverse  a  priori  covariance 
matrix  chosen  as  the  weighting  matrix  for  the  total  parameter  correction 
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term  in  the  cost  function.  A  minimum-variance  estimate  without  a  priori 
statistics  is  called  an  absolutely  unbiased  estimate;  when  a  priori  statistics 
are  included,  it  is  referred  to  as  an  unbiased  estimate. 

It  may  be  that  constraints  among  the  parameters  are  part  of  the  physical 
problem.  For  example,  suppose  that  precise  knowledge  of  the  relative  loca¬ 
tions  of  two  nearby  observing  stations  is  available.  If  such  locations  were 
among  the  parameters  m  a  differential  correction,  it  would  be  important  to 
constrain  the  parameter  corrections  to  preserve  the  relative  locations.  In 
TRACE66,  this  can  be  accomplished  by  introducing  linear  constraints  of  the 
form  P  =  BP7  +  C.  In  terms  of  corrections,  this  reduces  to  AP  =  BAP7  . 

The  weighted  least  squares  problem  may  be  solved  in  terms  of  AP7 ,  and  the 
above  matrix  equation  used  to  compute  AP. 

5.2  SEQUENTIAL  DIFFERENTIAL  CORRECTION  BY 

WEIGHTED  LEAST  SQUARES 

Systematic  repeated  application  of  batch  weighted  least  squares  is  generally 
referred  to  as  sequential  weighted  least  squares.  This  method  can  be  used 
with  or  without  including  a  priori  statistics.  However,  the  usual  manner  of 
processing  sequential  least  squares  includes  a  priori  statistics,  and  the 
result  is  interpreted  as  a  minimum-variance  estimate. 

In  the  computational  procedure,  the  data  is  separated  into  a  series  of  batches, 
or  stages.  After  each  stage  is  processed  iteratively  to  convergence,  the 
epoch  time,  vehicle  position  and  velocity  vectors,  any  remaining  parameters, 
and  the  parameter  error  variance/covariance  matrix  produced  by  the  estima¬ 
tion  process  are  updated  (propagated  forward)  in  time.  For  example,  the 
previous  epoch  quantities  are  updated  to  the  end  of  the  particular  data  batch. 
These  quantities  then  comprise  the  initial  conditions  and  a  priori  statistics 
for  the  next  data  stage.  In  this  fashion,  an  essentially  unlimited  amount  of 
data  over  an  unlimited  time  period  can  be  processed,  maintaining  dependence 
through  the  propagated  covariance  matrix. 
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The  a  priori  covariance  matrix  is  used  as  the  weighting  matrix  for  the 
parameter  correction  term  in  the  cost  function.  The  use  of  the  inverse 
normal  matrix  as  a  covariance  matrix  is  based  on  several  assumptions  which 
are,  in  general,  not  rigorously  satisfied.  Therefore,  an  optimistic  estimate 
of  the  parameter  error  variance/ covariance  matrix  is  produced,  and  syste¬ 
matic  errors  may  be  introduced  into  the  algorithm.  To  reduce  the  effects 
of  such  errors,  TRACE66  is  capable  of  using  additive  and  multiplicative 
deweighting  matrices  to  alter  the  a  priori  covariance  statistics  before  proces¬ 
sing  any  data  in  a  particular  stage.  This  is  normally  used  as  a  limited- 
memory  function,  so  that  older  data  measurements  have  less  influence  on  the 
current  solution. 

5.  3  POTTER  SQUARE- ROOT  FILTER 

The  weighted  least  squares  algorithms  are  particularly  valuable  when  the 
number  of  measurements  being  processed  in  one  update  interval  is  large  in 
comparison  to  the  number  of  parameters.  When  the  opposite  is  true,  it  may 
be  computationally  advantageous  to  use  a  recursive  estimation  algorithm  for 
an  orbit  determination  (Ref.  10).  For  this  purpose,  TRACE66  has  a  square- 
root  form  of  a  minimum-variance  filter  similar  to  that  originally  introduced 
by  Potter  (Ref.  11).  The  vehicle  trajectory  and  associated  partial  derivatives 
may  be  generated  analytically  (only  vehicle  position  and  velocity  may  be 
selected  as  parameters),  or  by  numerical  integration  (full  parameter  capa¬ 
bility).  In  addition,  the  program  may  be  instructed  to  periodically  switch 
between  the  Potter  square- root  filter  and  the  sequential  least  squares  algo¬ 
rithm  during  the  processing  of  a  given  data  set. 


6.  COVARIANCE  ANALYSIS  OVERVIEW 


The  performance  of  an  estimation  algorithm  is  best  evaluated  in  the  operational 
environment  for  which  it  is  intended.  Since  this  is  frequently  impractical, 
linearized  covariance  analyses  have  been  constructed.  This  essentially  entails 
a  determination  of  the  effects  that  specified  error  sources  have  on  the  precision 
of  the  estimation  procedure.  These  error  sources  may  be,  for  example,  inac¬ 
curacies  in  station  locations,  random  errors  in  measurements,  or  errors  in 
differential  equation  parameters.  Covariance  analysis  of  the  batch  and  sequen¬ 
tial  least  squares  estimators  is  available  in  TRACE66. 


Note  that  the  covariance  analysis  procedure  does  not  require  an  actual  orbit 
determination,  nor  does  it  require  actual  measurements.  It  does  require  the 
specification  or  simulation  of  the  tracking  system,  i.  e.  ,  tracker  locations 
and  characteristics,  visibility  constraints,  measurement  types  and  rates,  a 
list  of  the  parameters  to  be  studied,  and  a  reference  orbit.  A  linearized 
covariance  analysis  may  then  be  obtained  by  matrix  manipulations.  The  basic 
output  consists  of  covariance  matrices  as  a  function  of  time,  which  provide 
uncertainties  in  the  orbit  determination  (P)  parameters  and  in  the  satellite 
position  and  velocity  based  on  the  following: 


•  A  priori  information  on  the  P  parameters 

•  Uncertainties  in  the  measurements 

•  Uncertainties  in  a  set  of  nonestimated  (Q)  parameters 

6.  1  BATCH  WEIGHTED  LEAST  SQUARES 

COVARIANCE  ANALYSIS 


The  batch  weighted  least  squares  covariance  analysis  algorithm  is  summarized 
in  Fig.  10.  Note  that  because  weighted  least  squares  is  a  nonrecursive  pro¬ 
cedure,  the  P  parameter  covariance  matrices,  C(P)^  and  C(P)  provide 
uncertainties  in  the  P  parameters  at  a  fixed  time  only.  The  matrices  depend 
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•  FORM  OF  NORMAL  MATRIX:  A  WA  *  - ; - WHERE'- 

T  1  T 

L  Ag  WAp  •  Aq  WAq 


Ap  =  PARTIALS  OF 

MEASUREMENTS  WRT  P 

Aq  =  PARTIALS  OF 

MEASUREMENTS  WRT  Q 


•  VARIANCE/ COVARIANCE  MATRICES  OF  INTEREST: 


-  _ 

*  ATpWAp' 


ORBIT  DETERMINATION  PARAMETER  COVAR¬ 
IANCE  MATRIX  BASED  UPON  P- PARAMETER 
A  PRIORI  INFORMATION  C(P)0  AND  MEAS¬ 
UREMENT  COVARIANCE  MATRIX  2 


C{P)0+a  =  C(P)0+( —  )  C(Q)  P-PARAMETER  COVARIANCE  MATRIX  WHICH 

M  WQ/  \OUJ  INCLUDES  EFFECTS  OF  UNCERTAINTIES  IN 

0- PARAMETERS 

where:  (|)=.c(plp[A>q] 


NOTE  THAT: 

1.  VARIANCE/COVARIANCE  MATRICES  APPLY  TO  PARAMETER  VALUES  AT  EPOCH  TIME  ONLY 

2.  P-PARAMETER  A  PRIORI  INFORMATION  C(P)Q  IS  INCLUDED  IN  ATWA  BY  ADDING  C(P)"1  TO  [a*WA  ] 

3.  MEASUREMENT  UNCERTAINTIES  DEFINED  BY  2  ARE  INCLUDED  IN  ATWA  BY  SETTING  W^1 

4.  Q-PARAMETER  COVARIANCE  MATRIX  C(Q)  IS  PRESPECIFIED  AT  EPOCH 


•  VEHICLE  STATE  COVARIANCE  MATRICES  AS  A  FUNCTION  OF  TIME  ARE  OBTAINED 
BY  TRANSFORMING  THE  PARAMETER  COVARIANCE  MATRIX  AS  FOLLOWS 

CIXI  -m  CIP)P(^)T  ♦(£  ♦  fs  fp)  CIO,  ( *  ♦  is  |E)t 
Up  1  f\dpl  \dQ  dp  dQ  I  \  do  dp  dQ  I 


where: 


(  — )  8  (  — )  ARE  TRAJECTORY  PARTIAL  DERIVATIVES  (  SOLUTIONS  OF  THE 
\dPl  \dQ  I  VARIATIONAL  EQUATIONS) 


NOTE  THAT  THE  FIRST  TERM  OF  C  (X)  GIVES  THE  EFFECT  OF  MEASUREMENT  UNCERTAINTIES 
AND  THE  SECOND  TERM  GIVES  THE  EFFECT  OF  Q-PARAMETER  UNCERTAINTIES 

TO  TRANSFORM  C  (X)  TO  ANY  OTHER  COORDINATE  SYSTEM  R  APPLY  THE 
MATRIX  EQUATION 


‘'■"(SMS)' 


Fig.  10.  Orbital  Statistics  Via  Covariance  Analysis 
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on  the  amount  of  tracking  data  included  in  the  computation  and  will  therefore 
vary  with  time  even  though  they  relate  to  a  parameter  value  at  a  fixed  time. 

On  the  other  hand,  the  vehicle  state  covariance  matrices  C(X)  and  C(R)  give 
uncertainties  in  the  satellite  position  at  the  current  time. 

The  batch  least  squares  covariance  analyzer  can  be  operated  in  three  modes: 
real-time,  postflight,  and  update.  In  the  real-time  mode,  all  covariance 
matrices  are  based  on  tracking  data  up  to  and  including  the  time  of  the  output. 
Therefore,  the  matrices  are  representative  of  the  uncertainties  that  would 
result  from  processing  data  in  a  real-time  fashion.  In  the  postflight  mode, 
covariance  matrices  are  based  on  all  tracking  data  available.  Therefore, 
there  is  just  one  C(P)p  and  C(P)  +  .  The  satellite  state  vector  matrices 
represent  the  uncertainties  that  would  result  from  an  attempt  to  reconstruct 
a  trajectory  in  a  postflight  analysis  fashion.  Finally,  the  update  mode  consists 
of  updating  satellite  position  uncertainties  in  time  based  on  an  a  priori  covari¬ 
ance  matrix,  not  on  tracking  data.  Therefore,  the  satellite  state  vector  covari¬ 
ance  matrices  represent  the  uncertainties  in  satellite  positions  that  would  result 
from  a  simple  trajectory  generation. 

6.2  SEQUENTIAL  WEIGHTED  LEAST  SQUARES 

COVARIANCE  ANALYSIS 

A  sequential  least  squares  covariance  algorithm  takes  on  a  different  form 
because  of  the  periodic  updates  intrinsic  to  the  algorithm.  Two  matrices, 
corresponding  to  C(P)p  and  C(P)p+(^  must  be  maintained  and  propagated  accord¬ 
ing  to  the  updating  process.  The  notations  "C"  for  the  former,  and  "P"  for 
the  latter  are  used;  Q  =  C(Q),  the  Q  parameter  variance/covariance  matrix.  The 
lower-case  letters  "p"  and  nq:'  are  used  to  denote  the  parameter  vectors.  The 
algorithm  is  then  summarized  in  Fig.  11.  Note  the  multiplicative  and  additive 
deweighting  r-  r  •  ces  F  and  M  that  are  applied  to  C.  A  possible  use  of  this 
algorithm,  nof  applicable  to  batch  weighted  least  squares,  would  be  an  adaptive 
determination  of  the  best  F  and  M  for  a  particular  satellite  and  tracker  geometry. 
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INITIALIZE 


C  =  p  PARAMETER  CONVARIANCE  MATRIX  USED  BY  ALGORITHM 

P  =  ACTUAL  p  PARAMETER  CONVARIANCE  MATRIX 

S  =  0 

i  =  0 


Fig.  11.  Schematic  of  Sequential  Least  Squares 
Covariance  Analysis  Algorithm 
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7.  LUNAR  GRAVITY  FIELD  ANALYZER  OVERVIEW 


In  addition  to  the  TRACE66  trajectory  analysis  and  orbit  determination 
capabilities  previously  discussed  (see  Sec.  5),  a  particular  portion  of  the 
program  has  been  designed  and  developed  specifically  to  assist  in  the  dynam¬ 
ical  determination  of  lunar  gravitational  constants  and  in  the  reconstruction  of 
lunar  satellite  orbits  (Ref.  5).  The  ultimate  purpose  of  this  software  is  to 
obtain  a  lunar  gravity  model  from  the  JPL/NASA  Deep-Space  Network  two-way 
doppler  tracking  data.  The  force  model  can  be  expressed  in  terms  of  a 
spherical  harmonic  representation  of  the  central  body,  a  surface  distribution 
of  up  to  750  disks  or  point  masses,  and  gravitational  potential  due  to  the  earth 
and  sun. 

Within  TRACE66  proper,  the  key  features  of  this  option  are  as  follows: 

•  Mass  parameter  selection  algorithm  by  arc 

•  Printer  plot  of  doppler  residuals  versus  the  vehicle's 
selenographic  ground  trace 

•  Residual  and  data  edit  summaries 

•  Complete  orbit  determination  if  the  numer  of  parameters 
is  less  than  iOO 

•  Accumulation  of  a  large  (256  parameters/arc)  normal  matrix 
that  can  be  output  onto  an  auxiliary  file  with  specific  modifications 

Note  that  if  there  are  more  than  100  parameters  in  the  parameter  set  of 
interest,  the  merging  of  normal  matrices  and  the  solution  of  the  normal 
equations  is  performed  outside  of  TRACE66,  using  auxiliary  software  and 
techniques  written  explicitly  for  a  large,  sparse  linear  system.  This  software 
has  been  employed  along  with  the  TRACE66  normal  matrix  accumulation  to 
solve  a  system  of  the  approximate  order  2000. 
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